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ABSTRACT 

Utilizing inverse uncertainty quantification techniques, structural health 
monitoring can be integrated with damage progression models to form probabilistic 
predictions of a structure’s remaining useful life. However, damage evolution in 
realistic structures is physically complex. Accurately representing this behavior 
requires high-fidelity models which are typically computationally prohibitive. In the 
present work, a high-fidelity finite element model is represented by a surrogate model, 
reducing computation times. The new approach is used with damage diagnosis data to 
form a probabilistic prediction of remaining useful life for a test specimen under 
mixed-mode conditions. 


INTRODUCTION 

Structural Health Monitoring (SHM) is motivated by the idea that knowledge of an 
individual structure’s current health state increases reliability through a systematic 
reduction of uncertainty [1], However, SHM in this sense is inherently reactive, as 
damage must occur before diagnosis of the structural health provides any useful 
information. This can limit the utility of SHM, especially when applied to fatigue 
loading where remaining useful life (RUL) is of interest. A more comprehensive 
approach would be to utilize SHM systems to detect and quantify damage, but then 
extend this information through a prediction of how the damage will propagate and 
ultimately impact the structural RUL. This is the motivation for damage prognosis. 

In general, two primary challenges exist for damage prognosis. First, the problem 
is probabilistic in nature. SHM sensors cannot provide a deterministic assessment of 
the exact damage state, and the models used to predict how that damage progresses 
are, at best, approximations of the true underlying physics. Therefore, uncertainty 
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quantification has a central role in damage prognosis. Second, the evolution of damage 
occurs at multiple length scales [2,3]. However, the macroscale is of primary interest 
for SHM-aided prognosis, since detectable damage typically falls into this regime. 

There are many fracture mechanics-based, one-dimensional growth models for 
macroscale fatigue cracks [4]. Most of these models stem from the work of Paris et al. 
[5] in which the crack growth rate is a function of empirical parameters and either the 
stress intensity factor (SIF) or the strain energy release rate. For idealized geometries 
and boundary conditions, these values can be computed analytically. A great deal of 
work has been devoted to fatigue damage prognosis using these types of models for 
both composite and metallic stmctures; see [3,6] for examples. In the context of SHM, 
these works typically use sensors to acquire in-situ damage accumulation data that are 
then used to inversely quantify the uncertainty in model parameters via Markov Chain 
Monte Carlo (MCMC) techniques [7]. These uncertainties can then be propagated 
back through the models to form a probabilistic prediction of RUL [3,6,7] . 

For SHM to be an effective tool, however, real-world structures must be 
monitored. Damage evolution in these structures progresses in two or three 
dimensions with mixed-mode driving forces [8]. Finite element (FE) analysis is 
capable of capturing this behavior, but often at too high a computational cost for use 
with MCMC methods, which can require thousands to millions of model simulations 
to reliably sample from the target posterior distribution [7]. Parallel MCMC 
algorithms [9] constitute an active area of research that could eventually alleviate this 
burden. However, these methods depend on a set of Markov chains (i.e., each 
realization is dependent on its predecessor) and, therefore, still require an intractable 
number of simulations if utilizing FE analyses. 

Surrogate models have garnered a great deal of attention for their ability to 
quantitatively represent the primary features of high-fidelity models at a fraction of the 
computational cost [7,10]. Research into high-fidelity damage prognosis over the last 
five years has focused primarily on employing surrogate modeling techniques [10-14]. 
Sankararaman et al. [10,11] used the FE software ANSYS to develop surrogate 
models capable of returning an equivalent SIF for a given three-dimensional crack 
configuration and multi-axial loading condition based on a characteristic plane 
approach. Ling and Mahadevan then integrated this approach with SHM data to 
develop probabilistic fatigue damage prognostics for aluminum test specimens [12]. 
Hombal et al. expanded on this work by presenting a two-stage technique for planar 
approximation of a non-planar crack [13]. The semi-analytical approaches in [10-14] 
retain some advantages of high-fidelity modeling but simplify the fundamental growth 
mechanics and restrict potential crack shapes. 

More recently, Hombal et al. [14] offered an approach based on non-parametric 
representations of arbitrary crack fronts that were dimensionally reduced via Principal 
Component Analysis (PCA) [15] for use as input to a surrogate model. In this way, 
complex crack evolution can be modeled in a low-dimensional space, reducing 
computation times while retaining the mechanics of the non-planar growth. However, 
since the growth is internal to the surrogate model, the method is inherently dependent 
on a user-defined growth rate function. Without access to SIF information, this 
prohibits the consideration of the uncertainty in growth rate parameters, which are 
known to exhibit a high degree of scatter. 

In the present work, a flexible surrogate model for high-fidelity fatigue crack 
growth simulations is developed based on a separation of the dynamic and quasi-static 



aspects of damage evolution. The proposed method enables modeling arbitrary crack 
geometries and does not restrict access to growth parameters, making it well suited for 
probabilistic prognosis. The remaining sections of this paper are organized as follows. 
First, the prognosis method is presented over three subsections discussing the growth 
algorithm, the surrogate modeling approach, and the SHM-based inverse uncertainty 
quantification problem, respectively. Next, an experiment to validate the new method 
is developed, followed by the presentation and discussion of results. 


DAMAGE PROGNOSIS METHOD 
Fatigue Crack Growth Algorithm 

FRANC3D is a fracture mechanics code wherein arbitrary, three-dimensional, 
geometrically-explicit cracks are inserted into existing finite element meshes and 
grown via re-meshing. Stress intensity factors (SIFs) are calculated by evaluating the 
M-integral at the crack tip [16]. Unfortunately, a single crack growth simulation can 
take hours to complete, with the majority of this time consumed by the FE solution at 
each growth step. The goal of the presented work was to apply surrogate modeling 
techniques to develop an analog to FRANC3D which only takes seconds to run. 

To accomplish this goal, the FE solutions, and, thus, the need for a mesh, were 
replaced by a surrogate model. The FRANC3D approach to modeling crack growth 
could then be adapted to represent a crack parametrically as a collection of m points in 
space (Figure 1). In the present work, only the front was modeled, although, in a more 
general approach, the crack surface could be represented explicitly. Crack growth is 
calculated on a point-by-point basis. An initial crack front is supplied by the user, and 
then the surrogate model provides the mode I, mode n, and mode El SIFs at each front 
point. It is important to note that, by definition, these values correspond to a local set 
of orthogonal axes at each crack front point. These axes are defined by the crack front 
normal ((q) and tangent (e 3 ) unit vectors and their cross product ( e 2 ). 
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Figure 1: General illustration of the proposed fatigue crack growth process in an arbitrary volume. 



If subjected to mixed-mode conditions, the crack can change directions or “kink.” 
A ki nk angle, 9 kink , can be calculated via a user-defined algorithm. At each front 
point, 9 kink is used to rotate the local axes about e 3 , resulting in a set of m new growth 
directions. The SIFs are converted to an effective equivalent SIF range, A K ee , which, 
for the present work, was defined as K[ max - K[ min , where these two terms 
correspond to the projection of the SIFs perpendicular to the new growth direction at 
maximum and minimum load, respectively. Finally, crack front points are advanced 
along the growth directions by assuming a constant growth rate over a small number 
of cycles, AN, such that the crack growth magnitude Aa = da/dN x AN, where 
da/dN is the one dimensional growth rate. Therefore, the accuracy of the model is 
dependent on the number of cycles per growth step, which, without a mesh, can be 
sufficiently small. The crack growth rate is defined by a user-specified crack growth 
law. Herein, Walker’s modified Paris Law is used to incorporate TCratio effects [17], 


da A K ee ] 

dN ~ ° [(1 - /?) 1_m ' 


( 1 ) 


Here, C, m, and n are empirical parameters, and R = K[ min / K[ max . Any crack growth 
model can be used in place of (1). Once the crack growth increment and direction are 
defined, the current crack front can be propagated forward in time for AN cycles. The 
new crack front points are then compared to the geometric boundaries of the host 
component. A spline is fit through the front and trimmed anywhere the crack has met a 
free surface, as shown in Figure 1. A new set of m crack front points is then 
interpolated along the spline, comprising the new crack front. The crack growth 
process is iterated, and the crack front at each subsequent step is stored until a 
stopping condition is met (e.g., K, > K IiCritlcal ). 


Surrogate Model 


Surrogate modeling involves representing complex physical models as an input- 
output (10) relationship. In general, as the number of inputs to the surrogate model 
increases, the training set required to characterize the parameter space grows 
exponentially. The total number of parameters to describe a single crack front in three 
dimensions is 3m, resulting in an intractable 10 relationship. However, as in [14], 
Principal Component Analysis (PCA) can be used to solve this issue. PCA is a 
statistical technique for dimensional reduction of a dataset consisting of correlated 
variables where data are transformed to a new set of variables called the principal 
components. These components are uncorrelated and ordered such that the first few 
typically account for the majority of the variance in the original dataset [15]. PCA was 
used to reduce the dimensions of both the crack fronts and their associated SIF 
profiles. As an example, crack fronts with 3m =60 were dimensionally reduced to only 
two parameters while still accounting for 99% of the variance. 

Furthermore, reducing complexity in the 10 relationship itself is necessary for 
tractability. By isolating the crack growth mechanics and only using the surrogate 
model to capture the quasi-static evaluation of the SIFs, the modeled relationship 
becomes simple, direct, and well-defined. Crack geometry and loading serve as the 
inputs and the output is a PCA-reduced representation of the corresponding SIFs. The 
PCA reduction can then be inverted to recover the original SIF solution for a given 



geometry with a mean absolute percentage error on the order of 0.5% at each front 
point. Gaussian Process Regression (GPR) was chosen to build the surrogate model 
because of its ability to both capture this complex 10 behavior and to quantify the 
uncertainty in the surrogate model predictions [7]. For use with MCMC, the surrogate 
model requires training data that are dispersed over the infinite number of potential 
crack geometries, which is feasible due to three training-specific advantages of the 
proposed method. First, a crack growth simulation consisting of S steps produces S 
training data points since the 10 relationship is independent of cycle count. Second, 
training data can be generated in parallel, reducing the upfront computational burden. 
Finally, training data can be added if the initial training set is deemed insufficient. 

SHM and Inverse Uncertainty Quantification 


In SHM, damage diagnosis data collected over the life of an individual component 
or structure can be used to inform the prognosis through what is kn own as Bayes’ 
Theorem of Inverse Problems [7,18]. Assuming unbiased, independent and identically 
distributed measurement errors (i.e., £,~/V(0, a 2 )), the relationship between the model 
and the experimental observations is expressed as 


= fi(Q) + % i = l,...,k, (2) 

where Y t , fi(Q), and Q are random variables representing the experimental 
measurements, model response, and model parameters, respectively. The number of 
observations is denoted by k. The goal is to determine the posterior density of Q given 
the observed realizations, v obs , of Y t . Considering q to be the realizations of Q, Bayes’ 
Theorem can be used to formulate the inverse problem as follows: 


i ^obs) 
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where n(q, v obs ) is the posterior density of interest, n 0 (q) is the prior density in which 
any a priori knowledge of the parameters can be incorporated, and p is the number of 
parameters. The SHM data influence the problem through the likelihood, n(v obs \q). 
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The denominator of (3) is intractable for most problems, especially when p is large. 
However, Markov Chain Monte Carlo (MCMC) techniques, whose stationary 
distribution is the posterior density in (3), allow for a solution to the inverse problem. 


EXPERIMENTAL VALIDATION, RESULTS & DISCUSSION 

Two tension-tension fatigue crack growth experiments were conducted using 
edge-notched, A1 2024-T3 specimens. Adapted from work by Ingraffea et al. [19], 
holes were drilled (Figure 2) to induce localized mix-mode conditions, causing nearby 
cracks to kink. A constant amplitude stress of 5.95 ksi was applied with R = 0.1 at a 



frequency of 10 Hz. The first test was stopped once the crack growth rate exceeded 
1 x 1CT 5 in/cycle. This was assumed to correspond to a conservative end of life (EOF) 

i /o 

condition, K, > K I E0L , where K I E0L = 26.8 ksi-in “ was determined by recreating the 
final, as-measured crack in FRANC3D, applying the maximum load, and computing 
Kj. From the second experiment, a crack breaching the hole was determined to arrest 
for at least 3 million cycles before reinitiating, which was accounted for in RUF 
calculations. The observed fatigue crack path for the first specimen is presented in 
Figure 2. While no sensor-based SHM data were available for these experiments, 
visual crack tip measurements were taken over the specimen lifespan. To simulate the 
capabilities of guided wave interrogation, Gaussian white noise was added to the x and 
y measurements with variances, a , of 0.0004 and 0.002 in", respectively. The five 
measurements used for prognosis are represented by circled ‘x’ marks in Figure 2. 

To train the surrogate adequately over the wide range of potential crack starting 
locations, thirty initial y locations were simulated, uniformly spaced from 0.0 to 
0.3782 inches, each with a length of 0.08 inches. The origin is as shown in Figure 2. 
Training simulations were carried out using FRANC3D with median extension steps 
of 0.008 inches and a crack front template radius of 0.004 inches. In both the training 
and prognosis simulations, growth increments and kink angles were calculated using 
the Walker model (1) and Maximum Tangential Stress (MTS) theory [16], 
respectively. The final set of training data consisted of 3,000 first order relationships 
between a given crack front geometry and the corresponding SIFs. While the proposed 
methodology allows for the inclusion of loading as an input, this was not necessary in 
the present example. The training data were fit with a set of GPR models using the 
scikit-leam module for machine learning in Python [21]. Using the surrogate model, 
run times for crack growth simulations were reduced from approximately 3 hours to 
under 12 seconds - nearly three orders of magnitude faster. A verification of the 
surrogate model was conducted but will not be discussed further. 

The Bayesian inverse problem was solved utilizing the acquired SHM data and an 
adaptive MCMC algorithm based on the Python module PyMC [22]. A bum-in of 

5.000 samples was used to encourage sampling from the true posterior density, and 
then 10,000 samples were drawn. The resulting chains were thinned by retaining every 
tenth sample, and errors were assumed to be unbiased, independent and identically 
distributed. For illustrative purposes, only two model parameters were considered to 
be random variables: (1) the starting location of the crack, y ioc , and (2) the Walker 
exponential parameter, n. A global sensitivity analysis [23] would typically be 
conducted to inform these choices but was excluded for simplicity. Uninformative 
prior densities were assumed such that y ;oc ~l/(0. 01, 0.365) and n~f/( 1.0, 6.0). Updated 
parameter probability density functions (PDF) were inversely determined, from which 

1.000 model realizations were constructed using a Monte Carlo sampling method. 

These realizations were used to generate a normalized histogram for the RUF, as 

well as prediction and confidence bounds for the anticipated path, as seen in Figure 3a 
and Figure 3b, respectively. Since it is possible that the crack will grow into the hole 
and arrest, the histogram of predicted RUF in Figure 3a is bimodal (i.e., two distinct 
regions of probability exist). While the predicted mean RUF is less than the observed 
value by 77,537 cycles, the observed value still falls within the predicted PDF. This 
discrepancy is likely due to the noise level and sparse nature of the data. Furthermore, 
in Figure 3b, it appears that the predicted mean path agrees more with the data than the 
observed path, which was viewed as a validation of the proposed methodology. 
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Figure 2: Abaqus [20] FE model, experimental crack path, and SHM data, designated by the 'x' marks. 



Figure 3: (a) Histogram of RUL samples compared with the observed RUL and (b) the predicted crack 
path with 95% confidence and prediction bounds compared with the observed crack path. 


SUMMARY 

A new method for probabilistic prognosis of fatigue crack growth was 
demonstrated. Prohibitively expensive stress intensity factor computations were 
replaced by an efficient surrogate model trained via high-fidelity finite element 
simulations. The proposed approach was validated through fatigue crack growth 
experiments with induced, localized, mixed-mode conditions. Simulated SHM data 
were used to inversely quantify the uncertainty in model parameters, including those 
associated with the crack growth rate. These uncertainties were then propagated 
through the modeling framework to successfully predict the observed RUL and crack 
growth path. Model run times were reduced from around three hours to under 12 
seconds. Although the experiment was relatively simple, the modeling framework was 
developed generally, and the same techniques can be applied to a variety of more 
complex, three-dimensional SHM problems requiring high-fidelity prognosis. 
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